######################################################
#PA 2018 General Election Turnout Main Model.
# When Women Run, Voters Will Follow (Sometimes): 
##Examining the Mobilizing Effect of Female Candidates in the 2014 and 2018 Midterm Elections 
#By Safarpour, Wyckoff Gaynor, Rouse, and Swers #
######################################################

#clean the environment.
rm(list=ls()) 

#Packages.
library(mvtnorm) 
library(stargazer)
library(stats)
library(arm)
library(plyr)
library(Hmisc)
library(dplyr)

#setwd
setwd("/Users/ACS/Dropbox/When Women Run/Revision_PoliticalBehavior/R&R Part 2/Publication Docs/Replication Data and Code/")

#Step 1. Bring in fully recoded dataset.

PA<- read.csv("PA2018GeneralDataFinal.csv", header=TRUE)

#Step 2. Check length.

length(PA$V1)

#Step 3. Run Model.

m1<-glm(votedGeneral~ womangeneral+ female2+womangeneral*female2+
          GenerationZ+womangeneral*GenerationZ+
          Millennials+womangeneral*Millennials+
          GenerationX+womangeneral*GenerationX+
          SilentGeneration+womangeneral*SilentGeneration+
          rep + womangeneral*rep+
          otherparty+ womangeneral*otherparty+
          voted2016r+
          Proptrumpvote2016+totalreceipts, data=PA, family=binomial(link="logit"))
summary(m1)

stargazer(m1, out="PA2018GeneralCompetitivenessControl.html",
          type="html",
          covariate.labels = c("Female Democratic Candidate", "Female", "Generation Z",
                               "Millennials", "Generation X","Silent Generation", 
                               "Republican", "Other Party",
                               "Voted in 2016 General Election", "County Trump Vote 2016", 
                               "Total Money Raised",
                               "Female Democratic Candidate*Female",
                               "Female Democratic Candidate*Generation Z",
                               "Female Democratic Candidate*Millennials",
                               "Female Democratic Candidate*Generation X",
                               "Female Democratic Candidate*Silent Generation",
                               "Female Democratic Candidate*Republican",
                               "Female Democratic Candidate*Other Party"),
          dep.var.labels="Voted in PA 2018 General Election", column.labels = "",
          title = "Effects of candidate gender, generation, gender,  party, prior voting, county Trump vote share and total money raised on 2018 Midterm Turnout in PA",
          notes        = "Results from logistic regression. Standard errors in parentheses. Baseline age category: Baby Boomers. Baseline party: Democrats. Model excludes district 5, making the effect of female candidate only the effect for female Democratic candidates.", 
          notes.append = TRUE, notes.align = "l", digits=3, single.row = T)

sample <- PA[complete.cases(PA$votedGeneral, PA$womangeneral, PA$SilentGeneration, 
                            PA$Millennials, PA$GenerationZ, PA$GenerationX, PA$female2,
                            PA$rep, PA$otherparty, PA$voted2016r, PA$Proptrumpvote2016, PA$totalreceipts)==T,] # Regression sample
str(sample)

#Step 4. Compute effects.
preds <- predict(m1, type = "response")
options(digits=4)
summary(preds) #mean turnout.

n_draws <- 1000
set.seed(1714)
vcovm1<- vcov(m1)
coefm1<-coef(m1)
sim_coefs_m1 <- rmvnorm(n_draws, coefm1, vcovm1) #Specify the 1,000 simulated coefficients.
rbind(coef(m1), apply(sim_coefs_m1, 2, mean)) # Check they are close to original


ppwomen0m1s <- NULL
for (i in 1:n_draws) {
  ppwomen0m1s[i] <- mean(invlogit(sim_coefs_m1[i,1] + sim_coefs_m1[i,2]*0 + 
                                    sim_coefs_m1[i,3]*sample$female2 + 
                                    sim_coefs_m1[i,4]*sample$GenerationZ +
                                    sim_coefs_m1[i,5]*sample$Millennials +
                                    sim_coefs_m1[i,6]*sample$GenerationX +
                                    sim_coefs_m1[i,7]*sample$SilentGeneration+
                                    sim_coefs_m1[i,8]*sample$rep +
                                    sim_coefs_m1[i,9]*sample$otherparty+
                                    sim_coefs_m1[i,10]*sample$voted2016r+
                                    sim_coefs_m1[i,11]*sample$Proptrumpvote2016+
                                    sim_coefs_m1[i,12]*sample$totalreceipts+
                                    sim_coefs_m1[i,13]*sample$female2*0+
                                    sim_coefs_m1[i,14]*sample$GenerationZ*0 +
                                    sim_coefs_m1[i,15]*sample$Millennials*0 +
                                    sim_coefs_m1[i,16]*sample$GenerationX*0 +
                                    sim_coefs_m1[i,17]*sample$SilentGeneration*0+
                                    sim_coefs_m1[i,18]*sample$rep*0 +
                                    sim_coefs_m1[i,19]*sample$otherparty*0))
}
mean(ppwomen0m1s) # 

ppwomen1m1.s <- NULL
for (i in 1:n_draws) {
  ppwomen1m1.s[i] <- mean(invlogit(sim_coefs_m1[i,1] + sim_coefs_m1[i,2]*1 + 
                                     sim_coefs_m1[i,3]*sample$female2 + 
                                     sim_coefs_m1[i,4]*sample$GenerationZ +
                                     sim_coefs_m1[i,5]*sample$Millennials +
                                     sim_coefs_m1[i,6]*sample$GenerationX +
                                     sim_coefs_m1[i,7]*sample$SilentGeneration+
                                     sim_coefs_m1[i,8]*sample$rep +
                                     sim_coefs_m1[i,9]*sample$otherparty+
                                     sim_coefs_m1[i,10]*sample$voted2016r+
                                     sim_coefs_m1[i,11]*sample$Proptrumpvote2016+
                                     sim_coefs_m1[i,12]*sample$totalreceipts+
                                     sim_coefs_m1[i,13]*sample$female2*1+
                                     sim_coefs_m1[i,14]*sample$GenerationZ*1 +
                                     sim_coefs_m1[i,15]*sample$Millennials*1 +
                                     sim_coefs_m1[i,16]*sample$GenerationX*1 +
                                     sim_coefs_m1[i,17]*sample$SilentGeneration*1+
                                     sim_coefs_m1[i,18]*sample$rep*1 +
                                     sim_coefs_m1[i,19]*sample$otherparty*1))
}
mean(ppwomen1m1.s) #


effect.01m1 <-  ppwomen1m1.s-ppwomen0m1s
summary(effect.01m1) #  (0-1 women candidate effect)

# CIs and put results together in matrix
elementsm1 <- list(ppwomen0m1s, ppwomen1m1.s, effect.01m1) # I split these two lines so the code is clear.
lapply(elementsm1, summary) # See summary for each element of the list
cim1 <- do.call("rbind", (lapply(elementsm1, quantile, c(.025,.975)))) # Get quntiles for each

resultsm1 <- cbind(cim1[,1], lapply(elementsm1, mean), cim1[,2])
colnames(resultsm1) <- c("2.5", "Mean", "97.5")   
rownames(resultsm1) <- c("No Women Ran", "Female Candidate", "Effect No Women-Female Candidate")
resultsm1 #

#Drop what we no longer need.
rm(ppwomen0m1s, ppwomen1m1.s, effect.01m1, cim1)


# Store mean effect and 95% CI to visualize results of m1 later.
meaneffect_m1<-resultsm1[[6]]
upperci_m1<-resultsm1[[3]]
lowerci_m1<-resultsm1[[9]]


################Get Predictions by Generation M1###########

SG <- subset(sample, SilentGeneration==1)
BB <- subset(sample, BabyBoomers==1)
GX <- subset(sample, GenerationX==1) 
M <- subset(sample, Millennials==1) 
GZ <- subset(sample, GenerationZ==1)

w0.SGm1	<- NULL
w1.SGm1 	<- NULL
w0.BBm1	<- NULL
w1.BBm1 	<- NULL
w0.GXm1	<- NULL
w1.GXm1 	<- NULL
w0.Mm1	<- NULL
w1.Mm1 	<- NULL
w0.GZm1	<- NULL
w1.GZm1 	<- NULL


for(i in 1:n_draws){ 
  w0.SGm1[i] <- mean(invlogit(sim_coefs_m1[i,1] + sim_coefs_m1[i,2]*0 + 
                                sim_coefs_m1[i,3]*SG$female2 + 
                                sim_coefs_m1[i,4]*SG$GenerationZ +
                                sim_coefs_m1[i,5]*SG$Millennials +
                                sim_coefs_m1[i,6]*SG$GenerationX +
                                sim_coefs_m1[i,7]*SG$SilentGeneration+
                                sim_coefs_m1[i,8]*SG$rep +
                                sim_coefs_m1[i,9]*SG$otherparty+
                                sim_coefs_m1[i,10]*SG$voted2016r+
                                sim_coefs_m1[i,11]*SG$Proptrumpvote2016+
                                sim_coefs_m1[i,12]*SG$totalreceipts+
                                sim_coefs_m1[i,13]*SG$female2*0+
                                sim_coefs_m1[i,14]*SG$GenerationZ*0 +
                                sim_coefs_m1[i,15]*SG$Millennials*0 +
                                sim_coefs_m1[i,16]*SG$GenerationX*0 +
                                sim_coefs_m1[i,17]*SG$SilentGeneration*0+
                                sim_coefs_m1[i,18]*SG$rep*0 +
                                sim_coefs_m1[i,19]*SG$otherparty*0))
  w1.SGm1[i] <- mean(invlogit(sim_coefs_m1[i,1] + sim_coefs_m1[i,2]*1 + 
                                sim_coefs_m1[i,3]*SG$female2 + 
                                sim_coefs_m1[i,4]*SG$GenerationZ +
                                sim_coefs_m1[i,5]*SG$Millennials +
                                sim_coefs_m1[i,6]*SG$GenerationX +
                                sim_coefs_m1[i,7]*SG$SilentGeneration+
                                sim_coefs_m1[i,8]*SG$rep +
                                sim_coefs_m1[i,9]*SG$otherparty+
                                sim_coefs_m1[i,10]*SG$voted2016r+
                                sim_coefs_m1[i,11]*SG$Proptrumpvote2016+
                                sim_coefs_m1[i,12]*SG$totalreceipts+
                                sim_coefs_m1[i,13]*SG$female2*1+
                                sim_coefs_m1[i,14]*SG$GenerationZ*1 +
                                sim_coefs_m1[i,15]*SG$Millennials*1 +
                                sim_coefs_m1[i,16]*SG$GenerationX*1 +
                                sim_coefs_m1[i,17]*SG$SilentGeneration*1+
                                sim_coefs_m1[i,18]*SG$rep*1 +
                                sim_coefs_m1[i,19]*SG$otherparty*1))
  w0.BBm1[i] <- mean(invlogit(sim_coefs_m1[i,1] + sim_coefs_m1[i,2]*0 + 
                                sim_coefs_m1[i,3]*BB$female2 + 
                                sim_coefs_m1[i,4]*BB$GenerationZ +
                                sim_coefs_m1[i,5]*BB$Millennials +
                                sim_coefs_m1[i,6]*BB$GenerationX +
                                sim_coefs_m1[i,7]*BB$SilentGeneration+
                                sim_coefs_m1[i,8]*BB$rep +
                                sim_coefs_m1[i,9]*BB$otherparty+
                                sim_coefs_m1[i,10]*BB$voted2016r+
                                sim_coefs_m1[i,11]*BB$Proptrumpvote2016+
                                sim_coefs_m1[i,12]*BB$totalreceipts+
                                sim_coefs_m1[i,13]*BB$female2*0+
                                sim_coefs_m1[i,14]*BB$GenerationZ*0 +
                                sim_coefs_m1[i,15]*BB$Millennials*0 +
                                sim_coefs_m1[i,16]*BB$GenerationX*0 +
                                sim_coefs_m1[i,17]*BB$SilentGeneration*0+
                                sim_coefs_m1[i,18]*BB$rep*0 +
                                sim_coefs_m1[i,19]*BB$otherparty*0))
  w1.BBm1[i] <- mean(invlogit(sim_coefs_m1[i,1] + sim_coefs_m1[i,2]*1 + 
                                sim_coefs_m1[i,3]*BB$female2 + 
                                sim_coefs_m1[i,4]*BB$GenerationZ +
                                sim_coefs_m1[i,5]*BB$Millennials +
                                sim_coefs_m1[i,6]*BB$GenerationX +
                                sim_coefs_m1[i,7]*BB$SilentGeneration+
                                sim_coefs_m1[i,8]*BB$rep +
                                sim_coefs_m1[i,9]*BB$otherparty+
                                sim_coefs_m1[i,10]*BB$voted2016r+
                                sim_coefs_m1[i,11]*BB$Proptrumpvote2016+                                
                                sim_coefs_m1[i,12]*BB$totalreceipts+
                                sim_coefs_m1[i,13]*BB$female2*1+
                                sim_coefs_m1[i,14]*BB$GenerationZ*1 +
                                sim_coefs_m1[i,15]*BB$Millennials*1 +
                                sim_coefs_m1[i,16]*BB$GenerationX*1 +
                                sim_coefs_m1[i,17]*BB$SilentGeneration*1+
                                sim_coefs_m1[i,18]*BB$rep*1 +
                                sim_coefs_m1[i,19]*BB$otherparty*1))
  w0.GXm1[i] <- mean(invlogit(sim_coefs_m1[i,1] + sim_coefs_m1[i,2]*0 + 
                                sim_coefs_m1[i,3]*GX$female2 + 
                                sim_coefs_m1[i,4]*GX$GenerationZ +
                                sim_coefs_m1[i,5]*GX$Millennials +
                                sim_coefs_m1[i,6]*GX$GenerationX +
                                sim_coefs_m1[i,7]*GX$SilentGeneration+
                                sim_coefs_m1[i,8]*GX$rep +
                                sim_coefs_m1[i,9]*GX$otherparty+
                                sim_coefs_m1[i,10]*GX$voted2016r+
                                sim_coefs_m1[i,11]*GX$Proptrumpvote2016+
                                sim_coefs_m1[i,12]*GX$totalreceipts+
                                sim_coefs_m1[i,13]*GX$female2*0+
                                sim_coefs_m1[i,14]*GX$GenerationZ*0 +
                                sim_coefs_m1[i,15]*GX$Millennials*0 +
                                sim_coefs_m1[i,16]*GX$GenerationX*0 +
                                sim_coefs_m1[i,17]*GX$SilentGeneration*0+
                                sim_coefs_m1[i,18]*GX$rep*0 +
                                sim_coefs_m1[i,19]*GX$otherparty*0))
  w1.GXm1[i] <- mean(invlogit(sim_coefs_m1[i,1] + sim_coefs_m1[i,2]*1 + 
                                sim_coefs_m1[i,3]*GX$female2 + 
                                sim_coefs_m1[i,4]*GX$GenerationZ +
                                sim_coefs_m1[i,5]*GX$Millennials +
                                sim_coefs_m1[i,6]*GX$GenerationX +
                                sim_coefs_m1[i,7]*GX$SilentGeneration+
                                sim_coefs_m1[i,8]*GX$rep +
                                sim_coefs_m1[i,9]*GX$otherparty+
                                sim_coefs_m1[i,10]*GX$voted2016r+
                                sim_coefs_m1[i,11]*GX$Proptrumpvote2016+
                                sim_coefs_m1[i,12]*GX$totalreceipts+
                                sim_coefs_m1[i,13]*GX$female2*1+
                                sim_coefs_m1[i,14]*GX$GenerationZ*1 +
                                sim_coefs_m1[i,15]*GX$Millennials*1 +
                                sim_coefs_m1[i,16]*GX$GenerationX*1 +
                                sim_coefs_m1[i,17]*GX$SilentGeneration*1+
                                sim_coefs_m1[i,18]*GX$rep*1 +
                                sim_coefs_m1[i,19]*GX$otherparty*1))
  w0.Mm1[i] <- mean(invlogit(sim_coefs_m1[i,1] + sim_coefs_m1[i,2]*0 + 
                               sim_coefs_m1[i,3]*M$female2 + 
                               sim_coefs_m1[i,4]*M$GenerationZ +
                               sim_coefs_m1[i,5]*M$Millennials +
                               sim_coefs_m1[i,6]*M$GenerationX +
                               sim_coefs_m1[i,7]*M$SilentGeneration+
                               sim_coefs_m1[i,8]*M$rep +
                               sim_coefs_m1[i,9]*M$otherparty+
                               sim_coefs_m1[i,10]*M$voted2016r+
                               sim_coefs_m1[i,11]*M$Proptrumpvote2016+
                               sim_coefs_m1[i,12]*M$totalreceipts+
                               sim_coefs_m1[i,13]*M$female2*0+
                               sim_coefs_m1[i,14]*M$GenerationZ*0 +
                               sim_coefs_m1[i,15]*M$Millennials*0 +
                               sim_coefs_m1[i,16]*M$GenerationX*0 +
                               sim_coefs_m1[i,17]*M$SilentGeneration*0+
                               sim_coefs_m1[i,18]*M$rep*0 +
                               sim_coefs_m1[i,19]*M$otherparty*0))
  w1.Mm1[i] <- mean(invlogit(sim_coefs_m1[i,1] + sim_coefs_m1[i,2]*1 + 
                               sim_coefs_m1[i,3]*M$female2 + 
                               sim_coefs_m1[i,4]*M$GenerationZ +
                               sim_coefs_m1[i,5]*M$Millennials +
                               sim_coefs_m1[i,6]*M$GenerationX +
                               sim_coefs_m1[i,7]*M$SilentGeneration+
                               sim_coefs_m1[i,8]*M$rep +
                               sim_coefs_m1[i,9]*M$otherparty+
                               sim_coefs_m1[i,10]*M$voted2016r+
                               sim_coefs_m1[i,11]*M$Proptrumpvote2016+
                               sim_coefs_m1[i,12]*M$totalreceipts+
                               sim_coefs_m1[i,13]*M$female2*1+
                               sim_coefs_m1[i,14]*M$GenerationZ*1 +
                               sim_coefs_m1[i,15]*M$Millennials*1 +
                               sim_coefs_m1[i,16]*M$GenerationX*1 +
                               sim_coefs_m1[i,17]*M$SilentGeneration*1+
                               sim_coefs_m1[i,18]*M$rep*1 +
                               sim_coefs_m1[i,19]*M$otherparty*1))
  w0.GZm1[i] <- mean(invlogit(sim_coefs_m1[i,1] + sim_coefs_m1[i,2]*0 + 
                                sim_coefs_m1[i,3]*GZ$female2 + 
                                sim_coefs_m1[i,4]*GZ$GenerationZ +
                                sim_coefs_m1[i,5]*GZ$Millennials +
                                sim_coefs_m1[i,6]*GZ$GenerationX +
                                sim_coefs_m1[i,7]*GZ$SilentGeneration+
                                sim_coefs_m1[i,8]*GZ$rep +
                                sim_coefs_m1[i,9]*GZ$otherparty+
                                sim_coefs_m1[i,10]*GZ$voted2016r+
                                sim_coefs_m1[i,11]*GZ$Proptrumpvote2016+
                                sim_coefs_m1[i,12]*GZ$totalreceipts+
                                sim_coefs_m1[i,13]*GZ$female2*0+
                                sim_coefs_m1[i,14]*GZ$GenerationZ*0 +
                                sim_coefs_m1[i,15]*GZ$Millennials*0 +
                                sim_coefs_m1[i,16]*GZ$GenerationX*0 +
                                sim_coefs_m1[i,17]*GZ$SilentGeneration*0+
                                sim_coefs_m1[i,18]*GZ$rep*0 +
                                sim_coefs_m1[i,19]*GZ$otherparty*0))
  w1.GZm1[i] <- mean(invlogit(sim_coefs_m1[i,1] + sim_coefs_m1[i,2]*1 + 
                                sim_coefs_m1[i,3]*GZ$female2 + 
                                sim_coefs_m1[i,4]*GZ$GenerationZ +
                                sim_coefs_m1[i,5]*GZ$Millennials +
                                sim_coefs_m1[i,6]*GZ$GenerationX +
                                sim_coefs_m1[i,7]*GZ$SilentGeneration+
                                sim_coefs_m1[i,8]*GZ$rep +
                                sim_coefs_m1[i,9]*GZ$otherparty+
                                sim_coefs_m1[i,10]*GZ$voted2016r+
                                sim_coefs_m1[i,11]*GZ$Proptrumpvote2016+
                                sim_coefs_m1[i,12]*GZ$totalreceipts+
                                sim_coefs_m1[i,13]*GZ$female2*1+
                                sim_coefs_m1[i,14]*GZ$GenerationZ*1 +
                                sim_coefs_m1[i,15]*GZ$Millennials*1 +
                                sim_coefs_m1[i,16]*GZ$GenerationX*1 +
                                sim_coefs_m1[i,17]*GZ$SilentGeneration*1+
                                sim_coefs_m1[i,18]*GZ$rep*1 +
                                sim_coefs_m1[i,19]*GZ$otherparty*1))}
effect.SGm1	<- w1.SGm1 - w0.SGm1 
quantile(effect.SGm1, c(.025,.975)) #to get CIs; 
mean(effect.SGm1) #to get mean effect

effect.BBm1	<-  w1.BBm1 - w0.BBm1
quantile(effect.BBm1, c(.025,.975)) #to get CIs; 
mean(effect.BBm1) #to get mean effect

effect.GXm1	<-  w1.GXm1 - w0.GXm1
quantile(effect.GXm1, c(.025,.975)) #to get CIs;  
mean(effect.GXm1) #to get mean effect= 0.04254

effect.Mm1	<-  w1.Mm1 - w0.Mm1
quantile(effect.Mm1, c(.025,.975)) #to get CIs; 
mean(effect.Mm1) #to get mean effect= 

effect.GZm1	<-  w1.GZm1 - w0.GZm1
quantile(effect.GZm1, c(.025,.975)) #to get CIs; 
mean(effect.GZm1) #to get mean effect

#Drop what we no longer need.
rm(w1.GZm1, w0.GZm1, w1.Mm1, w0.Mm1,w1.GXm1,w0.GXm1, w1.BBm1 , w0.BBm1, w1.SGm1 , w0.SGm1, SG, BB, GX, GZ, M)

#Effect by gender subgroup.
w0.Menm1	<- NULL
w1.Menm1 	<- NULL
w0.Womenm1	<- NULL
w1.Womenm1 	<- NULL

Men <- subset(sample, female2==0)
Women <- subset(sample, female2==1) 

for(i in 1:n_draws){ 
  w0.Menm1[i] <- mean(invlogit(sim_coefs_m1[i,1] + sim_coefs_m1[i,2]*0 + 
                                 sim_coefs_m1[i,3]*Men$female2 + 
                                 sim_coefs_m1[i,4]*Men$GenerationZ +
                                 sim_coefs_m1[i,5]*Men$Millennials +
                                 sim_coefs_m1[i,6]*Men$GenerationX +
                                 sim_coefs_m1[i,7]*Men$SilentGeneration+
                                 sim_coefs_m1[i,8]*Men$rep +
                                 sim_coefs_m1[i,9]*Men$otherparty+
                                 sim_coefs_m1[i,10]*Men$voted2016r+
                                 sim_coefs_m1[i,11]*Men$Proptrumpvote2016+
                                 sim_coefs_m1[i,12]*Men$totalreceipts+
                                 sim_coefs_m1[i,13]*Men$female2*0+
                                 sim_coefs_m1[i,14]*Men$GenerationZ*0 +
                                 sim_coefs_m1[i,15]*Men$Millennials*0 +
                                 sim_coefs_m1[i,16]*Men$GenerationX*0 +
                                 sim_coefs_m1[i,17]*Men$SilentGeneration*0+
                                 sim_coefs_m1[i,18]*Men$rep*0 +
                                 sim_coefs_m1[i,19]*Men$otherparty*0))
  w1.Menm1[i] <- mean(invlogit(sim_coefs_m1[i,1] + sim_coefs_m1[i,2]*1 + 
                                 sim_coefs_m1[i,3]*Men$female2 + 
                                 sim_coefs_m1[i,4]*Men$GenerationZ +
                                 sim_coefs_m1[i,5]*Men$Millennials +
                                 sim_coefs_m1[i,6]*Men$GenerationX +
                                 sim_coefs_m1[i,7]*Men$SilentGeneration+
                                 sim_coefs_m1[i,8]*Men$rep +
                                 sim_coefs_m1[i,9]*Men$otherparty+
                                 sim_coefs_m1[i,10]*Men$voted2016r+
                                 sim_coefs_m1[i,11]*Men$Proptrumpvote2016+
                                 sim_coefs_m1[i,12]*Men$totalreceipts+
                                 sim_coefs_m1[i,13]*Men$female2*1+
                                 sim_coefs_m1[i,14]*Men$GenerationZ*1 +
                                 sim_coefs_m1[i,15]*Men$Millennials*1 +
                                 sim_coefs_m1[i,16]*Men$GenerationX*1 +
                                 sim_coefs_m1[i,17]*Men$SilentGeneration*1+
                                 sim_coefs_m1[i,18]*Men$rep*1 +
                                 sim_coefs_m1[i,19]*Men$otherparty*1))
  w0.Womenm1[i] <- mean(invlogit(sim_coefs_m1[i,1] + sim_coefs_m1[i,2]*0 + 
                                   sim_coefs_m1[i,3]*Women$female2 + 
                                   sim_coefs_m1[i,4]*Women$GenerationZ +
                                   sim_coefs_m1[i,5]*Women$Millennials +
                                   sim_coefs_m1[i,6]*Women$GenerationX +
                                   sim_coefs_m1[i,7]*Women$SilentGeneration+
                                   sim_coefs_m1[i,8]*Women$rep +
                                   sim_coefs_m1[i,9]*Women$otherparty+
                                   sim_coefs_m1[i,10]*Women$voted2016r+
                                   sim_coefs_m1[i,11]*Women$Proptrumpvote2016+
                                   sim_coefs_m1[i,12]*Women$totalreceipts+
                                   sim_coefs_m1[i,13]*Women$female2*0+
                                   sim_coefs_m1[i,14]*Women$GenerationZ*0 +
                                   sim_coefs_m1[i,15]*Women$Millennials*0 +
                                   sim_coefs_m1[i,16]*Women$GenerationX*0 +
                                   sim_coefs_m1[i,17]*Women$SilentGeneration*0+
                                   sim_coefs_m1[i,18]*Women$rep*0 +
                                   sim_coefs_m1[i,19]*Women$otherparty*0))
  w1.Womenm1[i] <- mean(invlogit(sim_coefs_m1[i,1] + sim_coefs_m1[i,2]*1 + 
                                   sim_coefs_m1[i,3]*Women$female2 + 
                                   sim_coefs_m1[i,4]*Women$GenerationZ +
                                   sim_coefs_m1[i,5]*Women$Millennials +
                                   sim_coefs_m1[i,6]*Women$GenerationX +
                                   sim_coefs_m1[i,7]*Women$SilentGeneration+
                                   sim_coefs_m1[i,8]*Women$rep +
                                   sim_coefs_m1[i,9]*Women$otherparty+
                                   sim_coefs_m1[i,10]*Women$voted2016r+
                                   sim_coefs_m1[i,11]*Women$Proptrumpvote2016+
                                   sim_coefs_m1[i,12]*Women$totalreceipts+
                                   sim_coefs_m1[i,13]*Women$female2*1+
                                   sim_coefs_m1[i,14]*Women$GenerationZ*1 +
                                   sim_coefs_m1[i,15]*Women$Millennials*1 +
                                   sim_coefs_m1[i,16]*Women$GenerationX*1 +
                                   sim_coefs_m1[i,17]*Women$SilentGeneration*1+
                                   sim_coefs_m1[i,18]*Women$rep*1 +
                                   sim_coefs_m1[i,19]*Women$otherparty*1))
}
effect.Menm1	<- w1.Menm1 - w0.Menm1 
quantile(effect.Menm1, c(.025,.975)) #to get CIs; 
mean(effect.Menm1) #to get mean effect.

effect.Womenm1	<-  w1.Womenm1 - w0.Womenm1
quantile(effect.Womenm1, c(.025,.975)) #to get CIs; 
mean(effect.Womenm1) #

#Drop what we no longer need.
rm(w1.Womenm1, w0.Womenm1, w1.Menm1, w0.Menm1, Men, Women )



#calculate and plot results by party.
DEM <- subset(sample, dem==1) 
OP <- subset(sample, otherparty==1) 
REP <- subset(sample, rep==1)

w0.DEM<-NULL
w1.DEM<- NULL
w0.OP<-NULL
w1.OP<-NULL
w0.REP<-NULL
w1.REP<-NULL

for(i in 1:n_draws){ 
  w0.DEM[i] <- mean(invlogit(sim_coefs_m1[i,1] + sim_coefs_m1[i,2]*0 + 
                               sim_coefs_m1[i,3]*DEM$female2 + 
                               sim_coefs_m1[i,4]*DEM$GenerationZ +
                               sim_coefs_m1[i,5]*DEM$Millennials +
                               sim_coefs_m1[i,6]*DEM$GenerationX +
                               sim_coefs_m1[i,7]*DEM$SilentGeneration+
                               sim_coefs_m1[i,8]*DEM$rep +
                               sim_coefs_m1[i,9]*DEM$otherparty+
                               sim_coefs_m1[i,10]*DEM$voted2016r+
                               sim_coefs_m1[i,11]*DEM$Proptrumpvote2016+
                               sim_coefs_m1[i,12]*DEM$totalreceipts+
                               sim_coefs_m1[i,13]*DEM$female2*0+
                               sim_coefs_m1[i,14]*DEM$GenerationZ*0 +
                               sim_coefs_m1[i,15]*DEM$Millennials*0 +
                               sim_coefs_m1[i,16]*DEM$GenerationX*0 +
                               sim_coefs_m1[i,17]*DEM$SilentGeneration*0+
                               sim_coefs_m1[i,18]*DEM$rep*0 +
                               sim_coefs_m1[i,19]*DEM$otherparty*0))
  w1.DEM[i] <- mean(invlogit(sim_coefs_m1[i,1] + sim_coefs_m1[i,2]*1 + 
                               sim_coefs_m1[i,3]*DEM$female2 + 
                               sim_coefs_m1[i,4]*DEM$GenerationZ +
                               sim_coefs_m1[i,5]*DEM$Millennials +
                               sim_coefs_m1[i,6]*DEM$GenerationX +
                               sim_coefs_m1[i,7]*DEM$SilentGeneration+
                               sim_coefs_m1[i,8]*DEM$rep +
                               sim_coefs_m1[i,9]*DEM$otherparty+
                               sim_coefs_m1[i,10]*DEM$voted2016r+
                               sim_coefs_m1[i,11]*DEM$Proptrumpvote2016+
                               sim_coefs_m1[i,12]*DEM$totalreceipts+
                               sim_coefs_m1[i,13]*DEM$female2*1+
                               sim_coefs_m1[i,14]*DEM$GenerationZ*1 +
                               sim_coefs_m1[i,15]*DEM$Millennials*1 +
                               sim_coefs_m1[i,16]*DEM$GenerationX*1 +
                               sim_coefs_m1[i,17]*DEM$SilentGeneration*1+
                               sim_coefs_m1[i,18]*DEM$rep*1 +
                               sim_coefs_m1[i,19]*DEM$otherparty*1))
  w0.OP[i] <- mean(invlogit(sim_coefs_m1[i,1] + sim_coefs_m1[i,2]*0 + 
                              sim_coefs_m1[i,3]*OP$female2 + 
                              sim_coefs_m1[i,4]*OP$GenerationZ +
                              sim_coefs_m1[i,5]*OP$Millennials +
                              sim_coefs_m1[i,6]*OP$GenerationX +
                              sim_coefs_m1[i,7]*OP$SilentGeneration+
                              sim_coefs_m1[i,8]*OP$rep +
                              sim_coefs_m1[i,9]*OP$otherparty+
                              sim_coefs_m1[i,10]*OP$voted2016r+
                              sim_coefs_m1[i,11]*OP$Proptrumpvote2016+
                              sim_coefs_m1[i,12]*OP$totalreceipts+
                              sim_coefs_m1[i,13]*OP$female2*0+
                              sim_coefs_m1[i,14]*OP$GenerationZ*0 +
                              sim_coefs_m1[i,15]*OP$Millennials*0 +
                              sim_coefs_m1[i,16]*OP$GenerationX*0 +
                              sim_coefs_m1[i,17]*OP$SilentGeneration*0+
                              sim_coefs_m1[i,18]*OP$rep*0 +
                              sim_coefs_m1[i,19]*OP$otherparty*0))
  w1.OP[i] <- mean(invlogit(sim_coefs_m1[i,1] + sim_coefs_m1[i,2]*1 + 
                              sim_coefs_m1[i,3]*OP$female2 + 
                              sim_coefs_m1[i,4]*OP$GenerationZ +
                              sim_coefs_m1[i,5]*OP$Millennials +
                              sim_coefs_m1[i,6]*OP$GenerationX +
                              sim_coefs_m1[i,7]*OP$SilentGeneration+
                              sim_coefs_m1[i,8]*OP$rep +
                              sim_coefs_m1[i,9]*OP$otherparty+
                              sim_coefs_m1[i,10]*OP$voted2016r+
                              sim_coefs_m1[i,11]*OP$Proptrumpvote2016+
                              sim_coefs_m1[i,12]*OP$totalreceipts+
                              sim_coefs_m1[i,13]*OP$female2*1+
                              sim_coefs_m1[i,14]*OP$GenerationZ*1 +
                              sim_coefs_m1[i,15]*OP$Millennials*1 +
                              sim_coefs_m1[i,16]*OP$GenerationX*1 +
                              sim_coefs_m1[i,17]*OP$SilentGeneration*1+
                              sim_coefs_m1[i,18]*OP$rep*1 +
                              sim_coefs_m1[i,19]*OP$otherparty*1))
  w0.REP[i] <- mean(invlogit(sim_coefs_m1[i,1] + sim_coefs_m1[i,2]*0 + 
                               sim_coefs_m1[i,3]*REP$female2 + 
                               sim_coefs_m1[i,4]*REP$GenerationZ +
                               sim_coefs_m1[i,5]*REP$Millennials +
                               sim_coefs_m1[i,6]*REP$GenerationX +
                               sim_coefs_m1[i,7]*REP$SilentGeneration+
                               sim_coefs_m1[i,8]*REP$rep +
                               sim_coefs_m1[i,9]*REP$otherparty+
                               sim_coefs_m1[i,10]*REP$voted2016r+
                               sim_coefs_m1[i,11]*REP$Proptrumpvote2016+
                               sim_coefs_m1[i,12]*REP$totalreceipts+
                               sim_coefs_m1[i,13]*REP$female2*0+
                               sim_coefs_m1[i,14]*REP$GenerationZ*0 +
                               sim_coefs_m1[i,15]*REP$Millennials*0 +
                               sim_coefs_m1[i,16]*REP$GenerationX*0 +
                               sim_coefs_m1[i,17]*REP$SilentGeneration*0+
                               sim_coefs_m1[i,18]*REP$rep*0 +
                               sim_coefs_m1[i,19]*REP$otherparty*0))
  w1.REP[i] <- mean(invlogit(sim_coefs_m1[i,1] + sim_coefs_m1[i,2]*1 + 
                               sim_coefs_m1[i,3]*REP$female2 + 
                               sim_coefs_m1[i,4]*REP$GenerationZ +
                               sim_coefs_m1[i,5]*REP$Millennials +
                               sim_coefs_m1[i,6]*REP$GenerationX +
                               sim_coefs_m1[i,7]*REP$SilentGeneration+
                               sim_coefs_m1[i,8]*REP$rep +
                               sim_coefs_m1[i,9]*REP$otherparty+
                               sim_coefs_m1[i,10]*REP$voted2016r+
                               sim_coefs_m1[i,11]*REP$Proptrumpvote2016+
                               sim_coefs_m1[i,12]*REP$totalreceipts+
                               sim_coefs_m1[i,13]*REP$female2*1+
                               sim_coefs_m1[i,14]*REP$GenerationZ*1 +
                               sim_coefs_m1[i,15]*REP$Millennials*1 +
                               sim_coefs_m1[i,16]*REP$GenerationX*1 +
                               sim_coefs_m1[i,17]*REP$SilentGeneration*1+
                               sim_coefs_m1[i,18]*REP$rep*1 +
                               sim_coefs_m1[i,19]*REP$otherparty*1))
}
effect.DEM	<- w1.DEM - w0.DEM 
quantile(effect.DEM, c(.025,.975)) #to get CIs; 
mean(effect.DEM) #to get mean effect

effect.OP	<-  w1.OP - w0.OP
quantile(effect.OP, c(.025,.975)) #to get CIs;   
mean(effect.OP) #to get mean effect

effect.REP	<-  w1.REP - w0.REP
quantile(effect.REP, c(.025,.975)) #to get CIs; 
mean(effect.REP) #to get mean effect

#Step 5.Store the Effects.
Election<-c(rep("2018 General",11))
Effect<-c(meaneffect_m1, mean(effect.Womenm1), mean(effect.Menm1), mean(effect.DEM), mean(effect.OP), mean(effect.REP), mean(effect.GZm1),mean(effect.Mm1),mean(effect.GXm1),mean(effect.BBm1),mean(effect.SGm1))
LowerCI<-c(lowerci_m1,  quantile(effect.Womenm1, c(.025)), quantile(effect.Menm1, c(.025)), quantile(effect.DEM, c(.025)),quantile(effect.OP, c(.025)),quantile(effect.REP, c(.025)),quantile(effect.GZm1, c(.025)),quantile(effect.Mm1, c(.025)),quantile(effect.GXm1, c(.025)),quantile(effect.BBm1, c(.025)),quantile(effect.SGm1, c(.025)))
UpperCI<-c(upperci_m1, quantile(effect.Womenm1, c(.975)), quantile(effect.Menm1, c(.975)), quantile(effect.DEM, c(.975)),quantile(effect.OP, c(.975)),quantile(effect.REP, c(.975)),quantile(effect.GZm1, c(.975)),quantile(effect.Mm1, c(.975)),quantile(effect.GXm1, c(.975)),quantile(effect.BBm1, c(.975)),quantile(effect.SGm1, c(.975)))
PAgeneral<-data.frame(Effect, UpperCI, LowerCI, Election) 
PAgeneral$`Registrant Subgroup`<-factor(c(1,2,3,4,5,6,7,8,9,10,11), levels=c(1,2,3,4,5,6,7,8,9,10,11),
                                        labels=c("Overall", "Women", "Men", "Democrats","Other\nParty", "Republicans","Gen. Z", "Millennials", "Gen. X", "Baby\nBoomers", "Silent Gen."), ordered = T)
#write effects to data frame.
write.csv(PAgeneral, "PA2018GeneralCompetitivenessControl.csv", row.names=F)
